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Abstract 

A two-dimensional model to solve compressible Navier-Stokes equations for the 
flow through stator and rotor blades of a turbine is developed. The flow domains 
for the stator and rotor blades, are coupled by the Chimera method that makes grid 
generation easy and enhances accuracy because the area of the grid that have high 
turning of grid lines or high skewness can be eliminated from the computational 
domain after the grids are generated. The results of flow computations show 
various important features of unsteady flows including the acoustic waves 
interacting with boundary layers, Karman vortex shedding from the trailing edge of 
the stator blades, pulsating incoming flow to a rotor blade from passing stator 
blades, and flow separation from both suction and pressure sides of the rotor 
blades. 

1. Introduction 

Most of computer codes currently used to simulate aerodynamics in compressors and turbines 
are based on the steady basis with a turbulence model However these time-averaged steady-state 
approaches are limited in their capabilities to compute some of the most critical features 
encountered in turbomachinery. Among their most serious shortcomings are their inability to 
predict phenomena associated with unsteady nature of the flows, and their inability to accurately 
simulate heat transfer. This is particularly evident when computing flow through a stator and rotor 
combination. Because of the non-uniform flow velocities from the stator exit plane, the moving 
rotor blades always encounter a time-varying, or pulsating, inlet flow. This constantly changing 
flow distribution of the inlet flow causes deviations from the steady inflow and changes the flow 
characteristics in the rotor blade passage. The importance of simulation of unsteady flows in the 
rotor blade passages has been recognized as a critical capability for performance analysis and 
design. 

The development of a code to analyze the three-dimensional transient flow through a stator and 
rotor is a difficult task for various reasons. To analyze the flow through a stator and rotor blades, 
it is necessary to include at least 5 to 6 blades in the computations because of differences in the 
blade pitches between the rotor and stator blades. The rotor blades are constantly moving relative 
to the stator blades. The computation takes significantly more cpu time not only because of the 
large number of blades to be considered but also because of the unsteady flow simulation. The 
development of more advanced codes requires both new mathematical methods for numerical 
computations and significantly more powerful computer resources. 

The present study [1] was initiated with the background thus stated. Two major objectives of 
the study are, first, (a) to examine some numerical methods suitable for the flow through stator 
and rotor rows, and second, (b) to analyze the effect of unsteady flow from the stator to the rotor. 
The present work is performed using a two-dimensional model. However, the experience and 
results of the investigation will be useful in developing three-dimensional models. 
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The remainder of this paper describes a computational model for the transient flow through 
stator and rotor blades based on the two-dimensional compressible Navier-Stokes equations and 
shows results of flow computations. We assume the pitch of stator blades is identical to that of the 
rotor blade. Therefore, only one stator blade that is stationary and one rotor that moves linearly to 
the negative direction of the y-coordinate are considered in the computations. For the direction 
tangential to the rotor movement, periodic boundary conditions are used for both rotor and stator 
flow channels. The body-fitted grids are used in both stator and rotor regions, but the grid for the 
rotor region moves with the rotor blade. The compressible Navier-Stokes equations are solved in 
each domain, and coupled to each other by the Chimera method[2]. 


2. Basic Equations and the Solution Methods 

The compressible Navier-Stokes equation for the stator domain is written as 
Qt "*■ E* + Fy = R x + Sy 

where 
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with 

^xx = [$• + 2|i.)u x + Xvy]/Re 
T X y = n(u y + v x )/Re 

Xyy = [(X + 2(l)Vy + Xu x ]/Re 
r = ut xx + vt xy + P(a 2 ) x 


( 1 ) 


( 2 ) 


(3) 


S = UT X y + VTyy + P(a 2 )y 

Here, p is the density, u and v are the velocity components, p is the pressure, e is the total energy, 
and all the variables are non-dimensionalized by 

x = x7L, y = y/L 
P = P/Po. u = G7a 0l v = v/ao 

P — P/(Po3o), e = e/(p 0 ao) / « 


H = p/Ho, X = \/X 0 
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t = aot/L, Re = p 0 aoL/^o 


Here the quantities with an overbar is the dimensional quantity, and those with subscript o is the 
value at the inlet. The variables x and y are coordinates, L is a reference length, (j. and X are 
viscosities, a is the sound speed The pressure is related to the conservative flow variables by 


p = (y-l)(e*^p(u 2 + v2)) 


(5) 


For the rotor domain that moves with the rotor blade, Eq.(l) holds but the variables are 
redefined by 
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*xx = [$• + 2|i)u x + Xv y ]/Re 

x X y = p(u y + Vx)/Re 

Xyy = [(X + 2|i)Vy + Xu x ]/Re 

S = UT X y + VTyy + P(3^)y 

The pressure is related to the conservative flow variables by 
p = (Y-l)(e-lp(u 2 + 7 2 )) 

Other nondimensional variables are 
x = xTL, y = y/L 

P = p/Po. u = v = v/ao 

P = p/(p<A>)> e = e/(p 0 ao) 


( 6 ) 


(7) 


( 8 ) 


M- = p/p 0 , X = XA 


and t and Re are common for the two domains. 


t = aot/L, Re = p 0 aoL/p o 


(9) 


At the interface between the stator region and the rotor region, the relations between the 
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variables are 


U = U, v = v- V, p = p, p = p 
where V is the nondimensional velocity of the rotor given by 
V=V/ao 

Introducing Eq.(4) and Eq.(8) into 
P = P 


yields 


e - ^p(u 2 + v 2 ) = e - IpdT 2 + v 2 ) 
or equivalently 

c = e -ip[v 2 -(v-V 2 )] 
or 

e =c+lp[(v + V) 2 -v 2 ] 


( 11 ) 


3. Computational Procedure 


The Navier-Stokes equations are transformed onto a rectangular computational domain, and 
solved by the implicit method based on the LU decomposition [4]. The Navier-Stokes code is 
based on the second-order accurate implicit method. However, the boundary condition along the 
trailing cut and the Chimera coupling both use the computed results from the previous times step. 
Therefore, the accuracy of the whole solution reduces to the first order in time at least locally. The 
gnd for each region is generated independently by the elliptic-hyperbolic method[3]. Two grid 
systems are used, one 38x220 points per blade, and another 55x658 points per blade. The grid 
pome along the trailing cut in the gird are not made continuous across the cut because otherwise 
very high skewness of grid occurs under the periodicity constraints of the grid. 

The inflow boundary conditions for the rotor domain are taken from the most updated 
computational results for the stator domain by the Chimera method. At the same time, the exit 
boundary conditions for the stator domain are taken from the most updated computational results 
from the rotor flow domain. The points along the inlet for the rotor regions are not same as the left 
S* . ro * or grid so the front area of the rotor grid with significant turning can be 
avoided. Likewise, the exist boundary for the stator domain does not coincide with the right 
boundary of the stator grid. 


The conventional method of coupling the two domains is as follows. The grid line which is 
vertical at the exit of the rotor domain would match the vertical grid line at the flow entrance of the 
rotor g^d. However, the rotor grid is moving with the rotor blades, so flow boundary conditions 
are exchanged between the two domains by means of one-dimensional interpolations. Since the 
rotor grid is constantly moving relative to the stator grid, the interpolation relations also constantly 
changes. With the conventional approach, the grid for the flow exit of the stator and the grid at the 
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entrance to the rotor should be carefully aligned. Whenever the distance between the stator and the 
rotor is altered, a new set of grids should be generated. Generally the grid lines in a c-grid for a 
turbomachinery blade have high turns at the top left and bottom left comer of the grid for a blade 
where the accuracy of the coordinate transformation is poor. 

With the Chimera method, a code developed for a single rectangular computational domain is 
applied for each of the flow domain belonging to one blade. Therefore, it is not necessary to match 
the grid line of the exit boundary of the stator domain with the grid line of the inlet boundary of the 
rotor domain. With the Chimera coupling, the front area of the rotor grid with high turning of grid 
lines may be removed from the actively used computational domain. Likewise the rear domain of 
the stator where high skewness of the grid tend to occur can be removed. The Chimera coupling 
requires two-dimensional interpolation to exchange boundary conditions between the stator and 
rotor. The interpolative exchanges of the boundary conditions disturb continuity of the field values 
-cross the boundaries of the two domains, but the conventional methods have the same problem. 
The Qumera method may be applied not only for the blade rows, but can be easily applied when 
additional flow domains are involved. This feature will become attractive in three-dimensional 
simulation particularly when the additional flow geometries such as inlet eye, volute, and exit 
diffuser are included (for a radial turbine). 


4. Illustration of Results 

The coarser grids for the stator and rotor are shown in Figures 1, where only every other grid 
lines that are parallel to the blade surface are omitted for clarity. The finer grids of 55x658/blade 
were generated by interpolation of the coarser grids (not shown here). 

Figures 2 illustrates density distribution at two different stator-rotor positions, which were 
developed as follows. The coarser grids of 38x220/blade were used for this flow computation. 
There are some overlapping between the computed flow domains for the stator and the rotor 
blades. The mid point in the overlapping area will be denoted by x=A. The density distribution is 
plotted for x<A, and the density distribution from the rotor calculation is plotted for x>A. 
Discontinuity of the density contours is due to the interpolating scheme to couple the two domains 
in the Chimera scheme. Note, however, that even if no smearing effect due to interpolation exists, 
the wake from the trailing edge of the stator blade docs not go into the rotor domain, because the 
plotting on the rotor domain is based on the relative coordinate. The grid size in the wake region of 
the stator is significandy smaller than that in front of the rotor blade. The difference of the grid 
sizes is also responsible to the discontinuity of the contours than Chimera coupling effects. 

An examination of Fig.2 reveals that shapes of density contours at these different time points, 
namely at different stator-rotor positions, are different. It is caused by the unsteady inlet flow to 
the rotor primarily because of changing relative position of the rotor to the stator blades. To show 
the effect of the moving relative positions, distribution of several flow quantities along a vertical 
line indicated in Figure 1 are plotted in Figure 3 in order of increasing time in cycle 8. The initial 
numerical transient of the solution such as propagation of acoustic waves disappears after two 
cycles. At the eighth cycle, the same time-dependent flow distribution repeats during each cycle. 

The computation with the finer grid of 55x658 was performed for the same geometry and 
physical constants as for the computation illustrated in Figure 2 and 3. Because of a massive 
amount of data from the finer grid computations, we have not yet been able to analyze the results to 
a full extent. However, snapshots of density distribution for a few selected time steps are 
illustrated in Figure 4. Figure 4a is a density distribution that occurs during the initial numerical 
transient and shows interaction between acoustic waves and boundary layers. Although the initial 
numerical transient is not a true physical phenomena, the flow behavior during that period still 
provides important information because similar transient should occur in case a sudden change of 
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inlet or exit flow conditions change or the load to the turbine impulsively changes. Figure 4b is for 
a later time than Fig. 4a but still in a late stage of the initial numerical transient. This figure shows 
Karman vortex shedding from the trailing edge of the stator blade, and flow separation from the 
suction side of the rotor blade. Figure 4c is for a time step after the rotor blades has passed stator 
blades seven times. The figure shows a flow separation at the lower surface of the rotor blade and 
a vortex being shed. A periodic shedding of vortices from the same location has been observed. 
At this time of writing it is not certain, however, if this vortex shedding phenomena is simply due 
to periodic change of the inlet flow to the rotor blade or induced by the Karman vortices from the 
trailing edge of the stator. But the latter is be highly likely. Karman vortices from the trailing edge 
are also observed in Figure 4c (Karman vortices are more clearly shown in the next figure). Figure 
4d (the number of contour lines are significantly increased to show details more clearly) shows 
Karman vortex street from the trailing edge of the stator blade at a different time step. Unlike the 
vortex shedding from a simple circular cylinder, the geometrical pattern of Karman Vortices from a 
turbine blade seem to be much irregular. 


5. Concluding Remarks 

Two dimensional modeling of the flow through stator and rotor rows based on the Chimera 
method has been described. The computational results show several important aspects of the 
unsteady flow through stator and rotor passages. The unsteady flow phenomena captured in the 
present computational simulation are not only interesting from academic curiosity but also worth 
further investigation of practical interest, particularly to understand the mechanism of heat transfer 
in turbomachinery. 
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Figure 1 Grids for the stator and rotor domains (only half of the points in 
the directions intersecting the blades and trailing cut are plotted) 
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cL Relative total energy 


igure 3 


Time-dependent variations of spatial distributions of field variables in front 
of the rotor blade leading edge (see A-A line in Figure 1) 





4a: dining the initial numerical transient 

Figure 4 Snapshots of density contours at different time steps (55x658 points/blade) 
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4b: a later stage of the initial numerical transient 
Figure 4 continued 
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4c: a time step long after the initial numerical transient 
Figure 4 continued 
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4 ± another time step long after the initial numerical transient 
(with a higher contour density than 4c) 

Figure 4 continued 
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4d\* repeat of 4d (with a much higher contour density) 
Figure 4 continued 
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